#include "atom.h"
#include "../m_tree/mTree.h"
#include "aa_table.h"
#include <string.h>
#include <stdlib.h>
#include <assert.h>

/** Take a pdb file on standard in, and create the following files:

    pdb_atoms
    pdb_atoms.fat
    pdb_atoms.mTree
    pdb_atoms.mTree.fat

    You'll probably want a perl script if you are doing this in bulk... 
    
    Also, if those files already exist, the result will be garbage or a core file 
    (if you're lucky.)

*/


int main () {

  char line[83];  // I'm liking the PDB all of a sudden. ;) 80 characters + \n + \0
  struct atom * a;
  int storeNumber;
  /*  struct mTree * tree = mTree_new_cached("pdb_atoms", "pdb_atoms.fat", 
					 "pdb_atoms.mTree", "pdb_atoms.mTree.fat", 
					 sizeof(struct atom), 99999,
					 (double(*)(void*,void*))&atom_distance);*/ // (The pdb is limited to 99999 atoms per file)
  struct store * oStore = store_open_cached("pdb_atoms", "pdb_atoms.fat", sizeof(struct atom), 500, 99999);
  struct mTree * tree   = mTree_open_cached(oStore, "pdb_atoms.mTree", "pdb_atoms.mTree.fat", 99999, (double(*)(const void*,const void*))&atom_distance);

  // To make everything line up, we need a placeholder in position 0.
 
  off_t new = store_new(oStore);
  a = store_read(oStore, new);
  strncpy(a->name, "STRT", 4);
  a->residue = -1;
  strncpy(a->element,"--", 2);
  store_write(oStore, new);
  store_release(oStore, new);

  while((fgets(line, 82, stdin) != NULL) && (strncmp("ENDMDL", line, 6))) { // Stop after the first model for now...

    if(!strncmp("ATOM", line, 4) ||   // ATOM field.
       !strncmp("HETATM", line, 6)    // HETATM field.
       ) {


      // if(line[21] == 'A' || line[21] ) {
	// Nuke all the pesky decimal points...
	int i;

	storeNumber = store_new(oStore);
	line[11] = '\0';
	assert(storeNumber == atoi(&line[6]));
	a = store_read(oStore, storeNumber);

	a->serialNumber = storeNumber;

	
	a->residue_type = aa_get_short_abbrev(&line[17]);

	a->altLoc = line[16];

	//	printf("%s", line);

	strncpy(a->name, &(line[12]), 4);

	a->chain = line[21];

	line[26]='\0';
	a->residue = atoi(&(line[22]));

	for(i = 34; i < 38; i++) { // X
	  line[i] = line[i+1];
	}
	line[38]='\0';
	a->x = atoi(&(line[30]));

	for(i = 42; i < 46; i++) { // Y
	  line[i] = line[i+1];
	}
	line[46]='\0';
	a->y = atoi(&(line[39]));

	for(i = 50; i < 54; i++) { // Z
	  line[i] = line[i+1];
	}
	line[54]='\0';
	a->z = atoi(&(line[47]));
	

	/*for(i = 57; i < 60; i++) { // Occupancy
	  line[i] = line[i+1];
	}
	line[60]='\0';*/
	
	
	for(i = 63; i < 66; i++) { // Temperature
	  line[i] = line[i+1];
	}
	line[66]='\0';
	a->temperature = atoi(&(line[61]));
	
	strncpy(a->element, &(line[76]), 2);

	//	atom_print(&a);
	
	store_write(oStore, storeNumber);

	mTree_insert(tree, storeNumber);
	
	store_release(oStore, storeNumber);

	a = 0;

	//    }
    } else if(!strncmp("TER", line, 3)) {  // Not ATOM or HETATM - If it's a TER, we need a placeholder.  Otherwise, ignore.
      new = store_new(oStore);
      a = store_read(oStore, new);
      a->serialNumber = new;
      strncpy(a->name, "TER ", 4);
      a->residue = -1;
      strncpy(a->element,"--", 2);
      store_write(oStore, new);
      store_release(oStore, new);

      a = 0;

    }
  }
  //  mTree_print(tree);
  //  printf("MTree size: %d\n", mTree_size(tree));

  mTree_close(tree);
  store_close(oStore);
  return 0;
}

/*
ATOM    982  H   LEU A 328     -14.230  -2.494   2.738  1.00  3.21           H 
*/
 
